Tunneling conductance in Superconductor/Ferromagnet junctions: a self consistent 

approach. 
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We evaluate the tunneling conductance of clean Ferromomagnet /Superconductor junctions via 
a fully self-consistent numerical solution of the microscopic Bogoliubov-DeGennes equations. We 
present results for a relevant range of values of the Fermi wavevector mismatch (FWM), the spin 
polarization, and the interfacial scattering strength. For nonzero spin polarization, the conductance 
curves vary nonmonotonically with FWM. The FWM dependence of the self-consistent results is 
stronger than that previously found in non-self-consistent calculations, since, in the self-consistent 
case, the effective scattering potential near the interface depends on the FWM. The dependence 
on interfacial scattering is monotonic. These results confirm that it is impossible to characterize 
both the the FWM and the interfacial scattering by a single effective parameter and that analysis 
of experimental data via the use of such one-parameter models is unreliable. 

PACS numbers: 74.50+r, 74.45. +c, 74.78.Fk, 72.25.Mk 
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I. INTRODUCTION 

Refinements in fabrication techniques that have oc- 
curred over the past twelve years have made it possi- 
ble to create devices that exploit and display the in- 
terplay between ferromagnetic and superconducting or- 
derings in superconductor/ferromagnet (SF) heterostruc- 
tures. These developments have raised the possibility of 
building devices that may be used to manipulate and de- 
tect spin polarized currents. Such devices, besides being 
of obvious considerable scientific interest, have technolog- 
ical applications for spintronicSfi It has been suggested 
that the introduction of an F layer into solid-state qubits 
may allow a stabilization of the state of the junction with- 
out the need for an external field.- Study of SF bilayer 
devices may illuminate the behavior of other systems that 
undergo spin/charge separation. 3 Spin/charge separation 
occurs when an electron is injected into a superconductor 
at the gap edge. The charge is absorbed by the conden- 
sate and quickly carried away by the supercurrent while 
its spin excitation remains^ The different spin orderings 
of the S and F layers lead to spin dependent transport 
effects, providing ways to study phenomena such as spin 
polarization, spin-diffusion lengths, and spin-relaxation 
timesi^ Other potential applications and device geome- 
tries are discussed elsewhere ^iLS. Earlier work is exten- 
sively reviewed in Ref. 0. The knowledge gained from SF 
spin probes can lead, if the results are properly analyzed, 
to a better understanding of spin transport properties in 
different materials and nanostructures. 

Thanks to the new technologies mentioned above, 
transport in clean SF bilayers is a very active subject 
of experimental research, but it remains a difficult topic 
theoretically. The problem of calculating the conduc- 
tance of an S/F bilayer has been solved^^ but only in a 
non-self-consistent manner, that is, by assuming the su- 
perconducting order parameter to behave as a sharp step 
function at the interface. This assumption neglects the 
proximity effect, the "leaking" of superconductivity into 



the ferromagnet (and, vice versa, that of the magnetism 
into the superconductor) that occurs in real systems. 
Thus, the proximity effect makes it impossible to identify 
a precise location where the superconducting correlations 
end and the ferromagnetic ones begin. Even with this 
simplification, the non-self-consistent approach is still far 
from trivial, as one must consider Andreev reflection^ 
the spin asymmetries due to the exchange interaction, 
and the spin coupling at the interface. Nevertheless, the 
stark fact remains that neglecting a self-consistent order 
parameter amounts to no less than omitting all influence 
of S/F proximity effects. We remedy this situation in 
this paper, where we show how the difficulties associ- 
ated with the proximity effect can be overcome and how 
a self-consistent calculation of the conductance can be 
achieved. The inclusion of the order parameter in a fully 
self-consistent manner allows us to fully account for the 
influence of the proximity effect, leading to correct, accu- 
rate and more interesting results at the worthwhile cost of 
having to solve a much more complicated problem. The 
self-consistent procedure results in the physical situation 
of a scattering potential that is dependent on the wave- 
functions of the scattered particles and holes. There is 
no analytic approach that will give a fully self-consistent 
pair amplitude for this problem: a numerical approach 
must be taken. The spatial variations in the pair ampli- 
tude complicate the details of Andreev reflection^ since 
the density of states (DOS) can no longer be described 
as having a well defined BCS gap^ and one must con- 
sider sub-gap states and even gapless superconductivity. 
Since there is no exact analytic solution for the spatial 
variations in the pair amplitude at an SF interface, there 
is no analytic solution for the DOS either, nor a fortiori 
for the conductance. In short, this is a complicated four- 
component scattering problem without an analytic form 
for the scattering potential. Despite these difficulties, the 
proximity effect must not be neglected, as it has, as we 
shall see, a strong influence on the results and on the way 
experimental data must be analyzed. 
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Early experimental worki^ was done on devices with 
an insulating oxide layer between the S and F layers, thus 
incorporating a superconductor/insulator/ferromagnet 
(SIF) tunnel junction. The stable oxide layer prevented 
the diffusion of ferromagnetic impurities into the super- 
conductor, making it possible to work in the clean limit. 
One disadvantage of SIF systems is that a tunnel barrier 
inhibits the proximity effect. Another disadvantage of 
the SIF junction is the reduced spin coherence lengths. 
While the electrons merely tunnel through the insulator, 
the strong binding to the lattice leads to spin decoher- 
ence. These tunnel junctions provide good information 
on the superconducting DOS, but it is dubious whether 
they provide a good picture of the state of the ferromag- 
net near the interface, and they will certainly not reveal 
the consequences of the proximity effect. Other work 14 
in a spin transistor geometry demonstrated the possibil- 
ity of growing SF interfaces over a large area with no 
insulating oxide layer. The purpose of that study was 
to explore transverse spin currents, not to examine bi- 
layer conductance at the SF interface. Also, it is not 
clear that there was no diffusion of magnetic atoms at 
the interface. More recen t 15 i 16 work has focused on bi- 
layer conductance in planar junctions with sharp inter- 
faces and clean materials. Characteristic peaks in the 
conductance were observed at bias voltages correspond- 
ing to Ao, the bulk superconducting gap energy, and an 
enhanced conductance at zero bias. 

Point contact bilayers offer a convenient method for 
studying relatively abrupt interfaces! 17 i 18 i 19 i 20 These de- 
vices are made by growing an oxide layer, thick enough 
to suppress tunneling, on top of a planar ferromagnet. A 
small hole is made through the oxide, and a superconduc- 
tor is grown on top. This geometry prevents diffusion of 
ferromagnetic atoms during growth and allows for a uni- 
form magnetic field. In this way, it is possible to study 
experimentall y 17 i 18 the conductance of clean SF bilayers 
in the ballistic limit. A quicker method for making point 
contacts mechanically places the tip of a sharpened su- 
perconducting wire on a bulk F sample. 6 While such a 
technique gives poor control over the size and shape of 
the contact, it is technologically a very desirable proce- 
dure as it leads to quickly obtained results, and therefore 
it may be used to probe the spin states of many differ- 
ent types of F materials. However, the analysis of the 
data, and of conductance spectroscopy experimental re- 
sults generally, has been hampered by incomplete under- 
standing of how to relate data and theory. We will ex- 
tensively address this issue in this work. We will be able 
to show that there is a clear influence of the proximity 
effect on transport at SF interfaces, which will enhance 
our understanding of these devices. This will also show 
that fully self-consistent studies are necessary to actually 
understand these point contact devices. 

In the case of strong interfacial scattering, Andreev re- 
flection and the proximity effect are suppressed. In the 
limit of a very large interfacial barrier, the conductance 
reflects essentially the superconducting DOS. If we work 



with small or medium interfacial barriers, as is the case 
here, Andreev reflection and the proximity effect become 
very important, leading to a much more interesting but 
also much more difficult problem. In the small barrier 
limit Andreev reflection can lead to an enhancement of 
sub-gap conductance^ since a single electron from F must 
excite a pair of electrons in S. This agrees with experi- 
mental observations of zero bias conductancei 15 i 16 ' 17 i 18 
The proximity effect changes the local DOS^i in the 
vicinity of the junction, introducing sub-gap states and 
even gapless superconductivity. 

Semiclassical methods, such as the Eilenberger or Us- 
adel equations^ i 22 ' 23 i 24 can give a reasonable approxima- 
tion to the conductance curves 2 ^ in dirty systems. How- 
ever, such models are not appropriate for clean systems, 
and can lead to spurious predictionsj 15 i 25 A common phe- 
nomenological approach is to define a current polariza- 
tion parameter— based upon the difference in spin-up 
and spin-down DOS in F. The portion of the transmis- 
sion coefficient due to the Andreev reflected hole (AR 
hole) is calculated for an SN interface. Since the AR 
hole must be in the opposite spin band of the incident 
electron, the AR hole coefficient is modified by a simple 
function of the polarization parameter. This is a reason- 
able phenomenological approach, but it lacks any micro- 
scopic justification, and is usually of limited success*^ 

There have been many attempts to derive a fully 
microscopic model to describe the conductance in SF 
systems These studies have used an 
abrupt approximation for the pair amplitude at the inter- 
face, thereby focusing on elastic scattering and Andreev 
reflection but neglecting the proximity effect. These 
studies do predict some correct qualitative features of 
the junctions. The author of Ref. remarks, and we 
agree, that only a proper consideration of the proxim- 
ity effect through self-consistent methods will give the 
correct quantitative features as well. A self-consistent 
study^ relying on a 2-D tight binding model has been 
performed but only to study spatial variations of the spin 
and charge currents parallel to interface. Therefore, there 
is a need for a fully self-consistent microscopic treatment 
of the conductance in SF bilayers to obtain quantitatively 
accurate results. 

In this work we show how to numerically calculate con- 
ductances of SF bilayers using a fully self-consistent solu- 
tion to the Bogoliubov deGennes (BdG) equations with 
a net current. We require a fully self-consistent solution 
to the BdG equations to properly treat clean inhomo- 
geneous systems in three dimensions £ We will analyze 
the resulting eigenfunctions, with proper boundary con- 
ditions, and generate transmission probabilities for par- 
ticles scattered from F to S via the Blonder-Tinkham- 
Klapwijk (BTK) method^ These transmission proba- 
bilities can then be used to calculate TV curves, from 
which we will calculate the conductances. We find good 
qualitative agreement with experimental results, includ- 
ing an enhanced conductance at zero bias. We show that 
the Fermi wavevector mismatch at the interface, the ex- 
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change field in the ferromagnet, and the interface barrier 
scattering all have a significant and independent effect 
on the shape of the conductance curves. This paper rep- 
resents the early fruits of a completely new technique, 
which we expect will eventually be used to study more 
complicated geometries and conditions. Our fully self- 
consistent solutions will allow for a very careful consid- 
eration of the influence of the proximity effect on the 
conductance, even for very thick F layers. 

In the next section, we describe in some detail our nu- 
merical methods and the procedures that we follow to 
extract the conductance as a function of applied volt- 
age. The results are presented and discussed in detail in 
Sect. III. A brief final section recapitulates our conclu- 
sions and points to future directions. 

II. METHODS 

The systems we study here are planar junctions made 
of SF bilayers with atomically smooth interfaces. We 



Where H = p 2 z /2m— Ep(z)+e± + U(z) is a single-particle 
Hamiltonian where p z /2m is the contribution to the ki- 
netic energy from motion in the z direction. The contin- 
uous variable ej^ is decoupled from the z direction but of 
course it affects the eigenvalues e„, which are measured 
from the chemical potential. We describe the magnetism 
by an exchange field h(z) which takes the value ho in the 
F material and vanishes in S. Within the superconducting 
layer, Ep{z) is equal to Eps, the Fermi energy of the S 
layers measured from the bottom of the band, while in the 
ferromagnet we have Ep{z) = Efm so that in the F re- 
gions the up and down band widths are Ep-\ = EpM + h$ 
and Epi = Epm — ho respectively. As pointed out in 
Ref. [^, one should not assume that Epu — Eps and 
we therefore introduce the dimensionless Fermi wavevec- 
tor mismatch parameter A = Epm / Eps = (kpM /^fs) 2 ■ 
We also introduce the dimensionless magnetic strength 
variable / by ho = Epj^I, The I — 1 limit corresponds 
to the "half-metallic" case. Interfacial scattering is de- 
scribed by the potential U(z) which we take to be of the 
form U{z) — H8{z — z ) where z is the location of the 
interface. The dimensionless parameter Hb = mH/kFM 
(everywhere in this paper h — 1) conveniently charac- 
terizes the strength of the interfacial scattering. The 
amplitudes u^(z) and v^z) can be written down from 
symmetry relations^ We will consider physically rele- 
vant values of A, /, and Hb- We measure all lengths in 
terms of the inverse of hps, the Fermi wavevector in the S 



assume that the layers are semi-infinite in the directions 
perpendicular to the interfaces (the x — y directions). We 
will take the bands to be parabolic, thus e± = k\J1ra, 
where k± is the wavevector in the transverse direction 
and £j_ is the energy corresponding to the x — y vari- 
ables. The superconductor is assumed to be an s— wave 
material. 



We use a numerical diagonalization of the micro- 
scopic Bogoliubov-deGennes^i equations for this inhomo- 
geneous system. Given a pair potential (order parame- 
ter) A(z), to be determined self consistently, the spin-up 
quasi-particle (u^zj) and spin-down quasi-hole (v^(z)) 
amplitudes obey the BdG equations: 



(2.1) 

I 

material. The dimensionless F width in the z-dimension 
is Dp = kpsdp, that of S is D$ = kpsds, and that of 
the entire sample is D. The dimensionless spatial coordi- 
nates are denoted by Z = kpsz. The pair potential A(Z) 
must be found through the self-consistent condition: 

MZ) = ^]T' (ul(Z)vi(Z) + ui(Z)vi(Z)) tanh . 

n 

(2.2) 

We will show that we can treat this problem as a plane 
wave scattering problem, with a scattering potential that 
exists over a finite region of space. The scattering poten- 
tial is, of course, obtained through self-consistent meth- 
ods. We will then look at the transmitted and reflected 
waves sufficiently far away from the interaction region 
that they, too, are plane w&vesM- In this case, the spa- 
tial extent of the scattering potential is governed by the 
proximity effect. Therefore, we must take the sample size 
large enough that the pair amplitude is zero over a large 
fraction of F, and it is approximately equal to its bulk 
value Ao over a large fraction of S. We call these regions 
in which the pair amplitude is approximately constant 
the asymptotic regions. We must also take the total sam- 
ple size to be very large so that the minimum wavevector 
in the problem, k m i n — allows us to approximate 
a continuum of incident plane waves. In principle, we 
should take the sample size to be infinite, but we must 
choose a finite value for computational considerations. 



H-h{z) A(z) 
A(z) -(H+h(z)) 



v [ n {z) 



v l n {z) 
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We have taken a total sample size of sixty times the su- 
perconducting coherence length £ > which turned out to 
be sufficient to avoid finite size effects. 

The procedure for calculating conductances begins by 
using the BdG equations (|2.f [) and the self consistency 
condition (|2.2j) to find a self-consistent order parameter 
for the system. We use a procedure similar to that in 
previous work4 i 12 i 36 i 37 We start with a fully three dimen- 
sional wavefunction 5'(r) = e t]i± ' r (u(Z),v(Z)) T . The 
factor of e lkj -' r contributes only e± = k\/2m, reducing 
this to a quasi-lD problem in the Z-direction. We then 
expand the u(Z) and v(Z) eigenfunctions in a basis of 
both sines and cosines: 



i±{Z) = 



'jyCOs(k q Z) 

% sm(k q Z) 



(2.3) 



where the ± signs in the left subindex refer to the sine 
or the cosine function respectively. This choice of basis 



is equivalent to using complex exponentials, but we gain 
some computational advantage in the very time consum- 
ing step of calculating a self-consistent pair amplitude by 
working with real rather than complex numbers for the 
time being. The wavevectors k q are defined in units of 
Ufs as: 



kq — 



2irq 

IT 



(2.4) 



with q being a positive integer. While the basis in prin- 
ciple requires all q > 0, in practice it is sufficient to chose 
a a cut-off large enough so that the largest wavevector 
in the problem corresponds to an energy that is a few 
ujd above the Fermi level. This choice of basis, which 
implies periodic boundary conditions as needed in this 
problem, allows for a wavefunction that is non-zero and 
has a non-zero first derivative at the boundaries, two nec- 
essary conditions for a non-zero current to be present. 



In this basis, the upper left quadrant 

- 



of the matrix in 

) pq + A(l-A-J) 



the left side of Eqns. (|2.ip is: 

sin((k p ~k q )D F ) _ sh\((k p + k q )D F ) 



k — k 
fvp rig 



and the lower right quadrant is: 

-2,o 



kp ~\-kq 



,D F ) + 1), 



-1)5, 



sin((fc p — k q )D F ) _ sin((k p +k q ) 
kp kg k p -\-kq 

1) 



for the cosine terms. 



;. For the sine terms: 



^B(cos(fe 9 D F ) 



TT+ = ± 

Q-P- JJ 



2 9 



) cos(k p Dp 

lU i 2 /i a j\ [ sin((k p -k q )D F ) sm((k p +k q ] 
) PQ' D V / k p — k q k p +k q 



and 



H, 



q-p- 



D y q 



- e± + l)5 pq - 



D V 1 



For the cross terms involving a sine 

_2_ / cos(( 



".p >* q 

-^M^( S m(k g D F ) S m(k p D F )), 

^ _|_ j\ sin((kp~k q )D F ) _ sin((k p +k q )D F ) 
' kp kq kp-\-kq 

^^{sm{k q D F )sm{k p D F )). 



u+ — 

U q+P- - JJ 



D 

and a cosine we have: 

-k q )D F ) _ cos((k p +k q )D)—cos((k p +k q )D F ) \ 
k p -\-k q J 



'k p -k q )D)-cos((k p 



t\ ( cos((k p ~k q )D F )-l _ cos((k p +k q )D F )-l \ 
V kp kq kp~\~kq J 



and 



-^M^{sm{k q D F )cos{k p D F )). 



H, 



q+p- 



2 / cos((fc p — k g )D)— cos((k p — k q )D F ) cos((k p -\-k q ) D) —cos((k p -\-k q ) D F ) \ 
^ k q kp~\~k q J 

's((fc p 



D 



D * ' \ kp k q 



cos((kp-\-kg)D F )-l ^ 
kp~\~kq 



\ n>p rLq 

+ WXHs [sin{k q D F ) cos(k p D F )). 



(2.5) 



(2.6) 



(2.7) 



(2. 



(2.9) 



(2-10) 
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In the above expressions, the matrix elements are given in u nits of Eps- i- e - m dimcnsionlcss form. 
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FIG. 1: (Color online) Scattering at an SF interface. An elec- 
tron is incident from the ferromagnet at an angle Q a with spin 
a and momentum k CT . Ordinary scattering leads to a partially 
reflected electron with the same spin and wavevector. For sub- 
gap scattering, excitations in the superconductor are Cooper 
pairs, consisting of an electron-like quasiparticle (ELQ) and a 
hole-like quasiparticle (HLQ), carrying a net spin of zero and 
charge 2e. To conserve charge, a hole is retro-reflected (An- 
dreev reflection) at an angle Q„ into the opposite spin band, 
denoted by a. 



The self consistency condition, Eq. (|2.2p , requires us 
to calculate the (<fip(Z)\A(Z)\(f> q (Z)) matrix elements 
numerically. The iterative procedure is then to cycle 
through Eqns. (|2.1[) and (|2.2[) until convergence is ob- 
tained, as has been explained elsewhere i 8 ' 12 i 36 ' 37 The very 
large sizes required for this problem, however, make the 
computations technically more demanding: it is neces- 
sary to massively vectorize all programs in order to ob- 
tain results in a reasonable amount of physical time. 

Once we have obtained a fully self-consistent spectrum 
of eigenfunctions and eigenvalues, we must extract a con- 
ductance from them. To successfully apply a plane wave 
scattering approach, we require a sample large enough so 
that we can assume that the scattering potential is con- 
fined to a finite region about the interface and that there 
are "asymptotic" regions in S and F, far away from the 
scattering region, where the scattering potential is not 
felt and the behavior is bulk-like. The non-self-consistent 
work 9 - assumes the scattering region to be infinitely small. 
The general outlines for the self-consistent and non-self- 
consistent treatments are similar. We will describe here 
first the salient features of the non-self-consistent treat- 
ment, and then discuss the changes necessary to study 
the self-consistent case. 

Consider a single electron of spin a in F with momen- 
tum k(j, which is at an angle 9 a with the Z-axis (see 
Fig. [1} . This electron may be partially reflected from the 
interface as an electron with spin a and a momentum of 



the same magnitude, but with the k z component in the 
opposite direction (ordinary reflection). Andreev reflec- 
tion allows a single charge from F to create a Cooper 
pair excitation in S without violating charge conserva- 
tion. In this process, when a single electron in spin band 
a is incident on the interface, a hole is retro-reflected into 
opposite spin band, which we denote by a, with momen- 
tum kg, and a Cooper pair excitation is created in S. We 
write this as: 



ik a Z 



+ a a e 



-ik„Z 





(2.11) 

where k a and k s indicate the if-components of k,j and k 5 , 
repectively. As noted above Eq. (|2.3p . k^ is conserved. 
We normalize to the incident particle flux, associate the 
ordinary reflected electron with the b a amplitude and 
the Andreev reflected hole with the a a amplitude. The 
transmitted wavefunction in S is: 



ik' 



Vs 



-ik'„z ( V<T 
U„ 



(2.12) 



where c a corresponds to an electron-like quasiparticle 
(ELQ) moving to the right and d„ corresponds to a hole- 
like quasiparticle (HLQ) moving to the left. The u a 
and v s amplitudes must obey the normalization condi- 
tion u a + v^ = 1 . Similar equations can be written down 
for states with incident holes. We can use the a a and b a 
amplitudes to write down a formula for the conductance 
in the T -> limits 



G(e,0) = G T M) + G X M) 



= TzZ p °( 1 + k r^ 2 -\ b ° 



,(2.13) 



where P a = (1 + p a I)/2 accounts for the probability for 
the incident electron to have spin with p a = 1 for the 
spin-up band and p a = — 1 for the spin-down band. The 
ratio y~ accounts for the different spin band wavevectors. 
The conductance is a function of the incident angle of the 
electron from F, 9 = 9 a . While we do not explicitly write 
it down, a a , b a , fc s , and k a are functions of e and 9 as 
well. The angularly averaged conductance is given by 9 



G a (e) = (G a ) - 



G a {e,9)cos(9)d9 
jj 1 ' cos{9)d9 



(2.14) 



where n, a is the angle of total reflection (critical angle) 
for incident particles of spin a . 

This model has the usual features of plane wave scat- 
tering. The conservation of k± across the interface leads 
to a modified 9 version of Snell's law: 



k a sin(0) = k s sin( 



(2.15) 
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and 



and 



fc ff sin(0) = fc;sin(0- 



(2.16) 



The incident and Andreev reflected angles are respec- 
tively given by sin(0) = k±/k a andsin(# 5 ) = k±/k&. The 
transmitted angle Q' a is found through sin(^) = k^_/k' a . 
This leads to a number of phenomena similar to those 
in electromagnetic wave scattering. In particular, the 
critical angle for ordinary reflection, which is given by 



a 



- 1 m 



= sin 



depends on the 



x /A(l+p„/) / 

spin band of the incident electron. There is also an angle 
beyond which Andreev reflection is no longer possible, 



given by £Iau — sin 



l + pal 



The foregoing discussion applies also to the non-sclf- 
consistent approach except that Eqs. (|2.1ip and (|2.12|) 
are now possible only in the asymptotic regions far from 
the interface. The comments that follow are the steps 
necessary to analyze the self-consistent results. Since 
there is no condition in the BdG equations to impose a 
net current traveling to the right the self-consistent spec- 
trum of eigenfunctions would not be, in the asymptotic 
regions, in the convenient form of Eqs. (|2.1ip and (|2.12|) 
even if we were to use a basis of complex exponentials. 
Our choice of periodic boundary conditions allows us to 
impose the condition of a net current a posteriori. This 
is very similar to the treatment given to one-dimensional 
scattering problems in elementary quantum mechanics. 
Since a one-dimensional problem with periodic boundary 
conditions will produce two- fold degenerate solutions, 38 
this quasi-one-dimensional problem also has two-fold de- 
generate solutions. We can find a linear combination of 
each pair of degenerate solutions that corresponds to a 
particle injected into the S layer from the F layer. If we 
look in the F layer sufficiently far away from the interface 
(in the asymptotic region) , so that the pair amplitude has 
gone completely to zero, we find that the numerically cal- 
culated eigenfunctions can be fit to the form: 



<P<r(Z) = 



r/u sin(fc (T Z + S a ) 
rj s sin(fc ff Z + 5 B ) 



Similarly, in the asymptotic region of S, we find 



(2.17) 



(2.18) 



If we take a degenerate pair of solutions, tpo-A(Z) and 
LPusiZ) , they will have the same k a and k s but differ- 
ent 77 CT , rjg-, 5 a , and 5 S so that they are orthogonal. We 
can then find a unique linear combination of a pair of 
degenerate states, ip a (Z) — Aip^AiZ) + BLp aB (Z), that 
is in the form of Eqs. (|2.11[) and (|2. 12[) in the asymptotic 
region. The appropriate complex A and B coefficients 
are 



A = 



it's A 



(2.19) 



B 



-Ms 



a A 



risBe 



: 



(2.20) 



where we have subscripted the parameters to indicate the 
eigenfunction to which they belong. A small amount of 
elementary algebra yields expressions for the a a and b a 
coefficients, which we can then substitute into Eq. (|2.13[) 
to find the conductance. 

To find quantities such as k a , k s , r) a , rjs, 5 a , and 5s, 
we analyze each eigenfunction in the asymptotic region of 
F. The wavenumbers, k a and kg-, are quantized in units 
of k m i n = 2-k/D, by the discretization of the system. 
To find their values in each case, we find the zeros of 
the eigenfunction to get the characteristic wavenumber 
of the eigenfunction and assign it to the nearest multiple 
of k m i n . It is then a trivial step to obtain values for S a , 
and 8s. Finally, we calculate the amplitudes, r) a and r/s, 
by integrating the square of the eigenfunction over an in- 
teger number of periods. We test the quality of our fit 
by integrating the fit equation with the numerical eigen- 
function over an integer number of periods. Each value of 
a a and b a depends on both the total energy value of the 
state, e, and the perpendicular wavevector k±. To ob- 
tain a smooth curve for the conductance, we must have 
several hundred values of e. Depending upon the critical 
angles, there will be tens to hundreds of k± corresponding 
to a single value of e. The above procedure must there- 
fore be repeated thousands of times to generate a single 
conductance curve. Since the a a and b„ coefficients, and 
therefore the conductance, are functions of these numer- 
ically obtained quantities, they are subject to numerical 
fluctuations arising from the discretization of the fc's and 
numerical errors. 

To reduce the influence of numerical fluctuations on 
the conductance, we use a procedure which has both nu- 
merical advantages and which better reflects the physics 
of the problem. At any temperature, the current is the 
integral over energies of the conductance: 39 



7(1/)= / G(e)(f(e-eV) - f(e))de 



(2.21) 



where G(e) — G-\(e) + Gj.(e) is the conductance, / is 
the Fermi function and V is the bias voltage. At low 
T, /(e) can be approximated by a step function, and it 
is immediately obvious that the elementary relationship 
between current and conductance 



G(V) 



dI(V) 
dV 



(2.22) 



is not just a definition but a mathematical identity. 

It is convenient to numerically evaluate G(V) as fol- 
lows: after the set of a's and b's have been obtained for 
all relevant energies and angles, we apply Eqs. (|2. 13|) . 
dUTU), and (j2T2Tj) to compute I(V). The integral in the 
first term of Eq. (|2.2ip (the second term is independent 
of V) is evaluated numerically by summing all computed 
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FIG. 2: (Color online) Example of the fitting procedure for 
the normalized current, as described in the text in connection 
with Eq. (|2.21fl . The example shown here corresponds to 
the case A = 1 and / = 0.866, which will be discused below 
(Figs.|4]and[5]). The crosses are the numerical results obtained 
as explained in the text and the curves are the fits. The units 
are explained in the text 



values of G(e) from an energy that is a few lob below the 
the Fermi energy up to the bias voltage, V. The exact 
value of the lower limit on the integral is not important 
because it will change the resulting value by a constant 
which will disappear when Eq. (|2.22p is applied. A non- 
zero lower limit is advantageous because the resulting 
value will depend on a larger number of states, thus re- 
ducing the effects of numerical fluctuations. With this 
procedure, the value of the integral depends on the con- 
ductance for all voltages below V. This has the effect of 
reducing the influence of numerical fluctuations, resulting 
in very smooth data, as illustrated by the example given 
in Fig. [21 The bias voltage is given throughout this work 
in units of Arj/e, so that a dimcnsionless bias voltage of 

V = 1 corresponds to the gap edge, and the normalized 
current follows from this and the dimensionless defini- 
tion of G discussed below. This resulting smooth data is 
then fit to an appropriate simple function. The fit for the 

V > 1 and V < I regions are done separately with the 
condition that the first derivative of the current be con- 
tinuous at V = 1. The resulting fit is then numerically 
differentiated to obtain the result for G(V). We follow 
this procedure for both the forward scattering and the 
angularly averaged cases. 



III. RESULTS AND ANALYSIS 

In this section, we present our numerical results for the 
conductance of the SF bilayers described above, at zero 
temperature. The materials that form the junction can 
be characterized by four dimcnsionless parameters: the 
exchange field parameter /, the Fermi wavevector mis- 
match A, the barrier strength Hb, and the dimensionless 
superconducting coherence length Sq = kps^o- We will 



vary the / and A parameters over their physically rel- 
evant ranges, which we take to be 0.1 < A < 2 and 
< I < 1 in our units. For the coherence length we 
choose So = 50 and consider both the case of negligi- 
ble barrier Hb = and that of an intermediate barrier, 
Hb = 1. For stronger interfacial scattering one of course 
quickly recovers the well-known and less interesting stan- 
dard tunneling results. Geometrically, we will use values 
of the dimensionless thicknesses Dp = D$ = 1500. By 
calculating conductance curves at a few values of A and 
/ for smaller and larger lengths, we have found that the 
lengths used here are sufficient to avoid finite size effects: 
using a smaller sample size increased numerical fluctua- 
tions although it did not change overall trends and aver- 
ages in the data. We therefore chose the largest sample 
size that allowed us to calculate high quality numerical 
results in a reasonable amount of computer time. 

Following common convention, dimensionless conduc- 
tances G(V) are normalized to e 2 /h. Thus, for a sample 
with A = 1, I = 0, Hb = and Ao = 0, (a homoge- 
neous non-magnetic conductor in the normal state) we 
would obtain G{V) = 1. The bias voltage is in units of 
A /e: thus V = 1 corresponds to the gap edge. We will 
present results for a range of values of A, both smaller 
and larger than unity, and consider values of the exchange 
field, I = 0.2, I = 0.5, and I = ^3/2 w 0.866. Results 
for G(V) will be given for two cases: forward scattering 
and angularly averaged scattering. For the forward scat- 
tering case, we integrate Eq. (|2.14p up to a small angle of 
7r/30 for spin- up and spin-down incident particles, while 
the angularly averaged results are integrated up to the 
maximum possible angle. Forward scattering results are 
most applicable to point contact devices. 

We consider first, in Fig. [31 the behavior of G(V) at a 
fixed intermediate value of I ( I = 0.2) for several values 
of A in the range of most experimental relevance, A < 1. 
This and subsequent figures are for Hb = 0, the experi- 
mentally more important case of a barrier free junction, 
which we will study in more detail. Later we will turn to 
the case where the tunneling limit is approached by set- 
ting Hb = 1. The forward scattering results are shown in 
the top panel of Fig. [3J Of the cases shown in this panel, 
those in the top curve (A = 1) serve as a good test of our 
technique, since even with a relatively small exchange 
field we find results qualitatively similar to the well es- 
tablished ones obtained by BTK22 at I = 0. In the BTK 
approximation, the conductance at zero bias voltage for 
samples with 1 = and A = 1 should be equal to 2 due 
to Andreev reflection, and we see that the top curve in 
our plot approaches this limit. There, the conductance 
at V — is slightly smaller than 2, which is easily ex- 
plained by a slightly suppressed Andreev reflection due to 
the different DOS in the spin-up and spin-down bands at 
our nonzero I. Furthermore, the conductance asymptot- 
ically approaches the normal state value (unity in these 
units) at larger bias voltages. We can therefore say that 
our results approach BTK for A = 1 and small I. 

The conductance is larger than unity throughout the 
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FIG. 3: (Color online) Results for the dimensionless conduc- 
tance G(V) (see text) as a function of bias voltage V, given 
in units of Ao / e. The top panel shows the forward scattering 
case for I = 0.2, and various A < 1. The forward conduc- 
tance is calculated by averaging over an incident particle flux 
with angles smaller than n/30. The bottom panel shows the 
angularly averaged case at the same values of I and A. The 
conductances are averaged over angles up to the critical angle, 
fie, for spin-up and spin-down particles. 



sub-gap range for all curves in Fig. [3] except at the small- 
est A. In every case, however, we see that the zero bias 
conductance is enhanced with respect to the large bias 
value G(V > 3). The value of G(0) monotonically de- 
creases with A, which is reasonable because a smaller A 
will lead to stronger ordinary scattering at the interface 
and inhibited Andreev reflection. If the Andreev reflec- 
tion and ordinary scattering responded to A in the same 
way, then the ratio of the zero bias conductance to the 
conductance at larger bias, say V = 3, would be con- 
stant. Instead, we see that the sub-gap enhancement is 
considerably reduced for smaller A, which implies that 
Andreev reflection is much more sensitive to A than or- 
dinary reflection. The ordering of the conductance curves 
for V > 1 is nonmonotonic. The order from greatest to 
least G is A = 0.71,1.0,0.50,0.25,0.11. This is sensible 



because the Fermi wavevector for the majority band in 
the F layer of the A = 0.71 curve is equal to 0.85 in units 
of the Fermi wavevector in S, while that for the A = 1.0 
curve is 1.2. Therefore, the A = 0.71 F layer is slightly 
better matched to the S than the A = 1.0 layer. 

The angularly averaged conductance (bottom panel of 
Fig. [3J is qualitatively similar to that found in the for- 
ward scattering case. The values of the conductances 
are smaller than the forward scattering case, because a 
smaller fraction of the incident current (to which every- 
thing is normalized) will pass across the junction. The 
zero bias enhancement (with respect to the large bias 
limit) is appreciably less pronounced at smaller A. At 
A = 0.11, the value of G(0) is not much larger than G(3). 
As in the forward scattering case, the values of G(l) de- 
pend strongly on A. The A = 1 curve in the bottom 
panel of Fig. [3] shows a very weak conductance peak at a 
sub-gap bias. This is a characteristic feature^ of systems 
for which the junction characteristics are dominated by 
the exchange field. For A = 0.71, which introduces a 
small amount of additional ordinary scattering, the con- 
ductance peak has completely disappeared. We shall see 
below that this peak is enhanced with larger I. 

The results for a very strong exchange field (I = 0.866) 
are shown in Figs. H] and [5j The first of these figures cov- 
ers the case where A < 1. The strong exchange field 
makes many features that are already present in Fig. [3] 
more dramatic. When present, the peaks at sub-gap 
biases are proportionally much higher. In the case of 
A = 0.25, there is now a very weak zero bias enhancement 
in the forward scattering case, and no such enhancement 
in the angular averaging case. For A = 0.11, the suppres- 
sion of the conductance for sub-gap bias voltages is very 
strong. As in the / = 0.2 case, the ordering of the con- 
ductance curves for V > 1 is nonmonotonic: the order 
from greatest to least is A = 0.50, 0.71, 1.0, 0.25, 0.11. As 
before, the ordering is dictated by the degree to which the 
Fermi wavevector of the majority band in F matches the 
Fermi wavevector in S. In all cases, the value is small be- 
cause the minority band contributes a minuscule amount 
to the conductivity. 

The next figure (Fig. [5]) is for / = 0.866, as in Fig. EJ 
but with A > 1. Such large values of A might be exper- 
imentally the case only if the superconducting material 
were some oxide material with s-wave pairing, but it is 
nevertheless of theoretical interest. For the two values 
of A > 1 considered there we see how the previously ob- 
served trends in A continue when the mismatch is in the 
opposite direction. We see that the subgap behavior is 
monotonic with a marked peak at nonzero bias, and that 
the conductance is nonmonotonic in A for V > 1. 

In Fig. we consider the effect of changing / at con- 
stant A. The mismatch is held to A = 1 and the ex- 
change field is taken from a moderate / = 0.2 to a strong 
I = 0.866. We chose to show the angularly averaged 
case because the sub-gap conductance peak is more ap- 
parent than in the forward scattering case. The I = 0.2 
curve, which was discussed above, appears nearly flat 
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FIG. 4: (Color online) Results for G(V) at a larger value of 
/ than in Fig. [3] / = \/3/2 ~ 0.866, and the same range of A. 
The top panel shows the forward scattering results and the 
bottom panel the angularly averaged values. The effects of a 
strong exchange field are readily apparent when this figure is 
compared with Fig. [3] (see text). 



when in the company of curves with larger I. For / = 0.5, 
the conductance peak is a little more pronounced. Over- 
all, the conductance is smaller as / increases because a 
larger exchange field leads to more poorly matched Fermi 
wavevectors from F to S. The I = 0.866 conductance 
curve is the most dramatic. The value at zero bias is 
proportionally much larger than G(l) than for / = 0.5 
or / = 0.2. The ratio of the zero bias conductance to the 
conductance at larger bias voltages is approximately two 
for all values of / shown here. This contradicts the com- 
mon assertion that the reduction in the minority band 
DOS invariably leads to a smaller Andreev reflection. 

In all figures shown thus far, the cusp in the conduc- 
tance for bias voltages corresponding to the gap edge 
is dependent on both A and /. The non-sclf-consistent 
analysis 9 predicts a value of G(l) that is dependent on / 
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FIG. 5: (Color online) Results for the same case presented in 
Fig. [4] and with the same panel arrangements, but for values 
of A in the range A > 1. See text for discussion. 
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FIG. 6: Angularly averaged values of G(V) for various /, at 
the same A = 1.0. This figure shows (see text) that stronger 
/ suppresses both ordinary scattering and Andreev reflection. 
Note that G(0)/G(3) ~ 2 for all values of I shown here, and 
that G(0)/G(1) increases with increasing I. 
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FIG. 7: (Color online) Plots of G(J; V = 1), the value of the 
dimensionless conductance at V = 1 in the forward scattering 
case for I = 0.2 (top panel) and I = 0.866 (bottom panel). 
The horizontal lines labeled by G(V, I) show the result from 
the non-self-consistent formula, Eq. (|3.1[) . The points rep- 
resent the numerical results from the present self-consistent 
study, at the values of I indicated in the legend. The error 
bars represent the uncertainty in the numerically obtained 
values (see text). Note the different horizontal range included 
in the upper and lower panels. 



but independent of A, namely: 



G(V = l;I) = 



4(1 — J 2 ) 1 / 2 
1 + (1-7 2 )V2 



(3-1) 



The non-self-consistent result depends very strongly on 
the assumption that the order parameter is independent 
of A, an assumption which is not valid, particularly near 
the interface, when the order parameter is calculated self- 
consistently. 37 Therefore one should look for violations 
of this relation in the correct, self-consistent calculation. 
Figure [7] is a direct comparison between the non-self- 
consistent formula 9 for G(V = I; I), Eq. (|3.1|) , and the 
self-consistent results. The figure shows G(V = 1;J) for 
I = 0.2 (top panel) and I = 0.866 (bottom panel). The 
results for / = 0.2 are plotted in the range of A included 
in Fig. [3] and those for I = 0.866 in the more extended 
range included in Figs. 2] and [5l The non-self-consistent 
results, which are, as explained above, indeppendent of 
A, are shown as the horizontal lines. The self-consistent 
results are the data points. The error bars represent the 
numerical uncertainty which arises because G(V) has a 
sharp cusp at V = 1. This discontinuous first deriva- 
tive at V — 1 makes this the most difficult value of the 
conductance to obtain numerically. Even with the nu- 
merical error, the difference between the two methods is 
apparent: there is a clear increasing trend in the self- 
consistent results as a function of A in contrast with the 
A-independent value predicted by the non-self-consistent 
result. This discrepancy is stronger at larger I and it af- 
fects (see bottom panel) both the small and large values 



of A, particualrly the latter ones. This large discrepancy 
here is due, as mentioned above, to the fact that the non- 
self-consistent result relies on the assumption that A is 
independent of A. 

In general, we find that for all values of I the exact 
results are in fair agreement with those obtained within 
the non-self-consistent approximation near A = 1 but 
that very considerable deviations occur in all cases for 
values of this parameter both smaller and larger than 
unity. Although in some intuitive way this may seem to 
make obvious sense, it is actually nontrivial, since the 
self-consistent result for A(z) near the interface is never 
equal to the non-self consistent (that is, the bulk) value. 
It seems, though, that the effect on this on the result 
for G(V = 1) is minimized when the mismatch is at a 
minimum. This is, however, rarely the experimental sit- 
uation: for most materials of experimental interest there 
is considerable mismatch. 

In the non-self-consistent treatment, it was found that 
the Fermi wavevector mismatch, A, has an important ef- 
fect on the conductance curves. 9 In the results found in 
this section, we have shown that the effect of A is even 
greater when self-consistently is included. It is a not un- 
usual experimental practic e) 40 ! 41 to characterize conduc- 
tance curves by using only a polarization parameter P, an 
interfacial barrier strength Hb, and Ao as fit parameters. 
In light of the results in Ref. |9j, confirmed and reinforced 
here, neglecting A is a deplorable practice, and can eas- 
ily lead to spurious results. A much sounder procedure 
would be to fit results for different samples using A, Hb 
and I as fitting parameters, taking Ao as the value of the 
bulk material and disentangling the effects of A and Hb 
by remembering that the latter, but not the former, will 
change from junction to junction. 

All of the previous results shown have been obtained 
in the regime of most experimental interest where the 
interfacial scattering parameter Hb is negligible. In the 
next figure, Fig. [8) we present results for Hb = 1, when 
the tunneling regime {Hb 3> 1) is approached, as can be 
seen from the very different shape of the curves. The two 
panels in the figure display both the dependence of the re- 
sults on / at constant A (top panel) and the dependence 
on A at constant I (bottom panel). In general, these 
dependences are not strong: this is because in the true 
tunneling regime G(V) simply reflects the superconduct- 
ing DOS and, as one approaches this regime, this should 
of course result in a weaker A dependence. This is well 
known to experimentalists: to obtain reliable measure- 
ments of quantities such as the polarization parameter 
it is always preferable to work with samples that have 
small interfacial scattering. In other words, a larger Hb 
inhibits Andreev reflection and the proximity effects. 

Nevertheless, at Hb = 1 some definite trends can be 
ascertained and definite statements can be made. In the 
top panel of Fig. [8] one can see that the I dependence 
at constant A is still relatively strong, particularly in the 
region V > 1. It is remarkable that the behavior with 
I is strongly nonmonotonic at larger voltages and also, 
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FIG. 8: (Color online) The conductance at interface scat- 
tering value Hb = 1, near the tunneling regime. All results 
correspond to the forward scattering case. The top panel 
shows G(V) for three values of I at A = 1, while the bottom 
panel shows G(V) at I = 0.866 for three values of A. 



although much more weakly, at V < 1. This makes in- 
terpolation schemes very doubtful in interpreting exper- 
imental data. In the bottom panel, corresponding to a 
strong magnet, we see that the mismatch dependence is 
weak in the region V > 1 for A < 1 but quite noticeable 
for A > 1. In the region V < 1 the situation is exactly the 
opposite: the curves corresponding to A = 1 and A = 2 
nearly coincide, while that for smaller A is clearly differ- 
ent. This nonmonotonic behavior contrasts again with 
that found in non-self-consistent results (see for example 
panel (b) of Fig. 2 in Ref. @) which, for similar values of 
J, vary monotonically 42 with mismatch in the same way 



over the whole range of V. 



IV. CONCLUSIONS 

We have introduced a method for calculating conduc- 
tances in an F/S bilayer within a fully self-consistent mi- 
croscopic model. Many of the features that we find are 
seen in experimental work: there is a sub-gap enhance- 
ment of the conductance, the conductance does approach 
the normal state value for larger bias voltages, and there 
is a cusp at bias voltages of unity. Most important, we 
show that, as already indicated by the non-self-consistent 
results, detailed experimental analysis (in particular the 
extraction of the spin polarization) is impossible if one 
does not take into account separately the effects of mis- 
match and those of barrier scattering. 

The features of the conductance curves agree only qual- 
itatively with those obtained via non-self-consistent pro- 
cedures. We find that there are strong quantitative differ- 
ences. One of the most obvious is that we find a strong 
dependence of G(V = Ao) on the parameter A, which 
characterizes the Fermi wavevector mismatch. In the 
non-self-consistent approach, G(V = Ao) is analytically 
found?, to be independent of A. In both treatments, the 
dependence of G{V = Ao) on I is similar in trend but 
quite different in detail. The conductance is reduced for 
larger A and larger I. We find a sub-gap conductance 
peak for strong I and A close to unity. All of this in- 
dicates that, while the non-self-consistent approach is a 
good tool to help us understand qualitatively some of the 
features of SF transport, a fully self-consistent approach 
is needed to properly model experimental data. 

This paper represents merely the first step in studying 
SF bilayers using a fully self-consistent pair amplitude. 
Future work may involve the addition of normal metal 
electrodes at the boundaries of the sample, allowing us 
to explore the effects of finite F and S widths. Detailed 
studies of the relationship between the local densities of 
states (DOS) and the conductance are also desirable. Fi- 
nally, the effects of finite temperature should be explored. 



Acknowledgments 

We are very grateful to Klaus Halterman for many dis- 
cussions on the technical aspects of this numerical work 
and to Igor Zutic for numerous conversations on this 
problem. 



* Electronic address: pbarsic@arete.com; Present address: 
Arete Associates, 1550 Crystal Dr. Ste. 703, Arlington, 
Virginia 22202 

t Electronic address: otvalls@umn.edu 



* Also at Minnesota Supercomputer Institute, University of 

Minnesota, Minneapolis, Minnesota 55455 
1 I. Zutic, J. Fabian, and S. Das Sarma, Rev. Mod. Phys. 

76, 323 (2004). 



12 



2 G. Blatter, V.B. Geshkenbein, and L.B. Ioffe, Phys. Rev. 
B 63 174511 (2001). 

3 S. Kashiwaya, Y. Tanaka, N. Yoshida, and M.R. Beasley, 
Phys. Rev. B 60, 3572 (1999). 

4 B. Leridon, J. Lesueur, and M. Aprili, Phys. Rev. B 72, 
180505(R) (2005). 

5 T. Yamashita, S. Takahashi, H. Imamura, and S. Maekawa, 
Phys. Rev. B 65, 172509 (2002). 

6 R.J. Soulen Jr., J.M. Byers, M.S. Osofsky, B. Nadgorny, 
T. Ambrose, S.F. Cheng, PR. Broussard, C.T. Tanaka, 
J. Nowak, J.S. Moodera, A. Barry, J.M.D. Coey, Science 
282, 85 (1998). 

7 A.I. Buzdin, Rev. Mod. Phys. 77, 935 (2005). 

8 K. Halterman and O.T. Vails, Phys. Rev. B 72, 060514 
(2005). 

9 I. Zutic and O.T. Vails, Phys. Rev. B 61, 1555 (2000). 

10 I. Zutic and O.T. Vails, Phys. Rev. B 60, 6320 (1999). 

11 A.F. Andreev, Sov. Phys. JETP 19, 1228 (1964). 

12 K. Halterman and O.T. Vails, Phys. Rev. B 69, 014517 
(2004). 

13 R. Meservey and P. Tedrow, Phys. Rep. 238, 173 (1994). 

14 M. Johnson, Appl. Phys. Lett. 65, 1460 (1994). 

15 S. Reymond, P. SanGiorgio, M.R. Beasley, J. Kim, T. Kim. 
and K. Char, Phys. Rev. B 73, 054505 (2006). 

16 S. Hacohen-Gourgy, B. Almog, and G. Deutscher, Appl 
Phys. Lett. 92, 152502 (2008). 

17 P. Raychaudhuri, A. P. Mackenzie, J.W. Reiner, and M.R. 
Beasley, Phys. Rev. B 67, 02041 1(R) (2003). 

18 S.K. Upadhyay, A. Palanisami, R.N. Louie, and R.A 
Buhrman, Phys. Rev. Lett. 81, 3247 (1998). 

19 F. Perez- Willard, J.C. Cuevas, C. Siirgers, P. Pfundstein 
J. Kopu, M. Eschrig, and H. v. Lohneysen, Phys. Rev. B 
69, 140502(R) (2004). 

20 P. Chalsani, S.K. Upadhyay, O. Ozatay, and R.A 
Buhrman, Phys. Rev. B 75 094417 (2007). 

21 K. Halterman and O.T. Vails, Physica C 397, 151 (2003) 

22 M. Krawiec, B.L. Gyorffy, and J.F. Annett, Phys. Rev. B 
66, 172505 (2002). 

23 M. Leadbeater, C.J. Lambert, K.E. Nagaev, R. Raimondi 
A.F. Volkov, Phys. Rev. B 59 12264 (1999). 

24 R. Seviour, C.J. Lambert, and A.F. Volkov, Phys. Rev. B 



59 6031 (1999). 

25 G.P. Pepe, R. Latempa, L. Parlato, A. Ruotolo, G. Au- 
sanio, G. Peluso, A. Barone, A. A. Golubov, Ya. V. Fomi- 
nov, and M.Yu. Kupriyanov, Phys. Rev. B 73, 054506 
(2006). 

26 M.J.M de Jong and C.J.W. Beenakker, Phys. Rev. Lett. 
74 1657 (1995). 

27 N. Stefanakis, Phys. Rev. B 64, 224502 (2001). 

28 R. Melin, Europhys. Lett. 51, 202 (2000). 

29 G. Tkachov, E. McCann, V.I. Fal'Ko, Phys. Rev. B 65, 
024519 (2002). 

30 J. Cayssol and G. Montambaux, Phys. Rev. B 71, 012507 
(2005). 

31 M. Bozovic and Z. Radovic, Phys. Rev. B 66, 134524 
(2002). 

32 K. Kuboki and H. Takahashi, Phys. Rev. B 70, 214524 
(2004). 

33 G.E. Blonder, M. Tinkham, and T.M. Klapwijk, Phys. 
Rev. B 25 4515 (1982). 

34 P.G. de Gennes, Superconductivity of Metals and Alloys. 
Addison- Wesley Publishing Co., Inc., 1966. 

35 E. Merzbacher, Quantum Mechanics, 2d ed, John Wiley & 
Sons, Inc., U.S.A.. 

36 PH. Barsic, K. Halterman, and O. T. Vails, Phys. Rev. B 
bf 75, 104502 (2007). 

37 K. Halterman and O.T. Vails, Phys. Rev. B 65, 014509 
(2002). 

38 J.W. Brown and R.V. Churchill, Fourier Series and 
Boundary Value Problems, 6th ed., McGraw-Hill, U.S.A. 
(2001). 

39 S. Datta, Electronic Transport in Mesoscopic Systems, 
Cambridge University Press, Great Britain (1995). 

40 L. Wang, TY. Chen, and C. Leighton, Phys. Rev. B 69, 
094412 (2004). 

41 L. Wang, TY. Chen, C.L. Chien, J.G. Checkelsky, J.C. 
Eckert, E.D. Dahlberg, K. Umemoto, R.M. Wentzcovitch, 
and C. Leighton, Phys. Rev. B 73, 144402 (2006). 

42 The mismatch parameter Lo used in Ref. |y is related to A 
as L\ = 1/A. 



